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FACILITY: MATH LIBRARY 
"ABSTRACT: 


MTHSDATANH returns the double precision hyperbolic arctangent of the 
double precision argument. The call is standard call-by-reference. 


VERSION: 2 
HISTORY: 
AUTHOR: 
Peter D Gilbert, 23-Jul-81: Version 2 
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‘ -SBTTL HISTORY ; Detailed Current Edit History 
3 ; VERSION 1 
8 1 $ 1-001 = Original from PL/I math Library 
009 § : Edit History for Version 02 of MTHSDATANH 
0000 5 : 2-000 Original 1981 
0000 6 ; 2-001 - = pone MOVZBL “y VTBL when accessing MTH$$AB_ALOG_V. PDG 2-Dec-1981 
0000 57 ; 2-00 hange RSB to RET after error exit. PDG 6-Jan-198 
4 28 ; 2-003 - od r problem with POLY instruction. PDG 19-Mar-1982 
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-SBTTL DECLARATIONS ; Declarative Part of Module 


4 11:21: 


MTHJACKET.MAR 


: EXTERNAL SYMBOLS: 


MTHSK 
MTHSSS1GNAL 
MTHS$AB_ALOG_V 
MTHSSAB_D_FHT 


GLOBAL 
BINARY,CALLS, CONDI TIONALS, DEF INITIONS EXPANSIONS 
c INVARGHAT 


: EQUATED SYMBOLS: 


OPDEF X, OP, SH 
ADDX *XO00@SH+0P,R'X,M'X 
ADDX3 *XO1Q@SH+0P,R'X,R'X,W'X 
SUBX *X02Q@SH+OP,R'X,M'X 
SUBKS “*XO3@SH+OP,R'X,R'X,W'X 
M *X04QSH+0P,R°X,M'X 
MULX3 *XOS5@SH+0P,R'X,R'X,W'X 
“X06@SH+0P,R'X,M'X 
DIVX3 “XO7@SH+OP,R'X,R'X,W'X 
CVTWX “XOD@SH+OP,RW,W'X 
POLYX “X15@SH+OP,R'X,RW,AB 
MOVX *x007D,RQ,WQ ; MOVQ 
MOVAX *X007E,AQ,WL ; MOVAQ 
D, <*x0060>, 0 
ADDX_  *X00a0+*x0060.RD.MD 
ADDX3 “*x01a senbhe »RD,RD,WD 
UBX nies +*X0 RDM 
SUBX3 “*x03a0+*xX $6 eRD,RD,WD 
MULX *x04a0+°X -RD,M 
MULX3 “x05a0+*xX0060,RD,RD,WD 
IVx *x0680+*X0060,RD.MD 
DIVX3 “x07a0+*x0060,RD,RD.WD 
CVTWX  “*xX0Da0+*xX0060,RW, 
POLYX %X15a0+*x0060,RD,RW,AB 
Ovx *x007D,RQ,WQ ; MOVQ 
MOVAX *X007E,AQ,WL ; MOVAQ 
7 ; Bit offset to exponent 
7 ; Bit offset to exponent 
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! 
0000 1 -PSECT _MTHSCODE PIC,SHR,LONG,EXE,NOWRT 

i 3 program section for math routines 
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| 
i 
| ~SBTTL MTHSDATANH - Double Precision Hyperbolic Arctangent 
| : FUNCTIONAL DESCRIPTION: 
; ; DATANH = Double precision floating point function 
; ; DATANH(X) is computed using the following approximation technique: 

If iX! >= 1.0, error. Otherwise 

Let (1#X)/(1-xX) = f * (2e*n), where 1/2 <= f < 1 


If n is gracter than or equal to 1 then 
1 =n- 1 and Fil = 2*f, 
se 


set N =n and fF = f. 


| 
| 
| 
If if e-li< qrens then 
2tatanh(X) = Neln(2) + Wo + WeP(W), 
mars ¢ = ((14F)/(1-F)) #2eeN = 1 
and P is. a polynomial of degree f= =5,D0=9. 


se obtain FHI arta C1 sqyet to F) from table loo 

; 2eatanh(X) = Ln (14#K)/(1-x)) = Neln(2) + entrniy” + 2*Q(2*2), 
; upere Bp s R 9 pehyogalat of of aneree F=2,D=5, 

; where Z = 


(Desen)#C19K) 7 OK) 


di variety of ways: 
1 + ¥"(1-X)J 


; 

X + Y = X*¥J 
y +X = X*¥] 
d 


*te 


+ X2(14¥))/ 
NOTE: The quantities Ln(A=FHI 
equations in two parts - a high part (containing the 
ht order al and a low part (containing the low 
order bits. In the code the a and low parts of the 
constants are indicated by a _HI] and _LO suffix i sonal 
tively. The values were chosen such Ehat N*LN2_HI + 
LN_FH “Hl is exactly representable. 
; ; CALLING SEQUENCE : 
datanh.wd.v = MTHSDATANH(x.rd.r) 
INPUT PARAMETERS: 
K=4 ; Contents of x is the argument 
: IMPLICIT INPUTS: none 


: OUTPUT PARAMETERS: 


+ X#(1-¥)) 


( 
1 
X + X*Y¥J/ | 
) ln2 are used in the above 


an 


00000004 
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MTHSDATANH ; Double Precision Hyperbo EP=1984 01:15: AX/VMS Macro V04-00 Pa 
5-008 MTHSDATANH = Double Pecteten Hyperbolic oS Ei obe 94339 9 EMTHRTL . SRERTHDATANH. MAR: 1 wa 
0 2 8 
1 : VALUE: Double precision hyperbolic arctangent of the argument 
0 18 ; IMPLICIT OUTPUTS: none 
0088 5 ; COMPLETION CODES: none 
0088 é ; SIDE EFFECTS: 
0088 24 > Signals: MTH$K_INVARGMAT if !X! >= 1.0 with reserved operand in RO (copied to 
0088 5 ; the signal mechanism vector CHFS$L_MCH_RO/R1 by LIBSSIGNAL). 
ppee § ; Associated message Wy "Floating Overflow in math Library'’’. Result is 
08 ; reserved operand -0.0 unless a user supplied (or any) error handler changes 
+ 8 ; CHFSL_MCH_RO/R1. 
0088 5 0 : NOTE: This procedure disables floating point underfiow and integer 
0088 1 ; overflow, causes no floating overflow or other arithmetic traps, and 
Bpes $36 ; preserves enables across the call. 
0088 234 ; Note: This routine is written to avoid causing any integer overflows, 
0088 235 ; floating overflows, or floating underflows or divide by 0 conditions, 
0088 236 ; whether enabled or not. 
0088 237 ; 
0088 238 ;--- 
0088 44 
0135 = 311 ones sey ERR: BRW ERROR 
O1FC 0088 seg -ENTRY MTHSDATANH, ACMASK ; standard call-by-reference entry 
008D 24 ; disable DV (and FU), enable IV 
50 04 BC 7D baat ste MOVX @aX(AP), RO ; RO = arg 
52 08 50 43 0091 ie SUBF3 RO, S*#1.0, R2 ; R2 = 1-X (okay to trunc towards 0) 
Fl 15 0095 247 BLEQ ERR 3 ATANH(X) is not defined for X>=1 
56 08 50 41 Q097 248 ADDF3 RO, S*#1.0, R6 ; R6 = 14+X (okay to truc towards 0) 
EB 15 O009B 249 BLEQ ERR ; ATANH(X) is not defined for X<=-1 
56 52 46 009D 250 DIVF2 R2, R6 ; R6 = approximation to (1+#X)/(1-x) 
58 00000000'GF 9E& QOAO 251 MOVAB G*MTHS$$AB_ALOG_V, R8 
58 68 CO OQ0A7 526 ADDL (R8), ; R8 = address of ALOG table 
55 56 OO7F 8F AB OQOAA = 25 BICW #iaF _EXP-1, R6, RS ; RS = Biased exponent 
55 4000 8F ry 00B0 254 SUBW 6. #*X4000, R5 : R5 = Unbiased exponent 
64 «1 Hey $22 BLEQ NEG_EXP _-¢ Branch to processing for n=<0 
55 0080 8F ag 0087 57 SUBW #1aF_EXP, R5 «Exponent is positive, RS =N=n- 1 
56 55 A geet 58 SUBW 5, R6 ; RO =F = 
56 56 9A OQOBF 59 MOVZBL R6, R6 ; R6 = index into ALOG table 
00000000 00C 60 IF NE, F_EXP=X_EXP : 
Bae 2] fs WIQ<F _LEXP=X_EXP>, R5 ; Shift R5 to scale X-floating 
7E 55 6D 00C § CVTWX RS, -(SP) ; Push N onto the stack 
58 6846 98 00C 64 CVTBL (R8)CR6], RB ; RB = offset into FHI tables 
5D 19 00C9 65 BLSS LN_1 PLUS_W ; Branch to handle F close to 1 
58 00000000'GF48 /7E bose 66 MOVAX G*ATASSAB_D_FHICR8], RB ; RB = Address of FHI 
56 88 7D 4 of MOVX (R8)+, R6 ; R6 = FHI 
006 $9 : Compute 2 = (F_= FHI)/(F + FHI) 
0D6 0; Z = Se - ig ptt Sa + ¥"(1-Xx) 
00D6 ‘i 3 Z=(1+X-Y¥ + Xev)/[1 +X + ¥ = X*Y 
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DATANH = Double Precision Hyperbolic 6-SEP=19 MTMRTL.SREIMTHDATANH.MAR: 1 Tab 
i if : where Y = FHI*2**N, roughly equal to (1+#X)/(1-x) 
56 55 Ad 00D : ADDW , Rb : R6 = FHI *® 2eeN = SFHI 
52. (08 33 4g 88 75 SUBX3 R6, $*#1.0, R2 + R2 = 1 = SFH 
56 (08 =~ «60 Dp 7 ADDX A#1.0, R6 + R6 = 1 + SFH HW] 
52 5666 E 7 DIVX. R6, R + R2 = (1-SFHI)/(1+SFHI) = D 
56 COS 50 61 OOE 78 ADDX3 RO, Re, R6 : R6=D + X 
5 5064 OE? 79 MULX RO : R2 =D * X 
52 «208 «=©60 «(O0EA 6 280 ADDX  $*#1.0, R2 + R2 = 1 + Dex 
56 (52)—Si«6 QED 81 DIVX R2, RO + R6 = (D+X)/(14D#xX) = 7 
G0r 388 ; Compute Z**2, P(Z**2) and Z2*P(Z**2) 
50 56 56 OFO 285 MULX3 R6, R6, RO 3 RO = Ze* 
FFSE CF 46005)6=6 500s 75s«OOF4 SBS POLYX RO, #LOGLEN2, LOGTAB2 : RO = P(Z**#2) 
50 56 OOFA 387 MULX _R6, RO : RO = Z2*P(Z**2) 
beep $80 ; Compute B = N*LN2_LO + LN_FHI_LO + Z*P(Z*2) 
52 FFO6 CF 6€ 65 OOFD $34 MULX3 (SP), LN2_LO, R2 ; R2 = N®*LN2_LO 
52 88 60 0103 332 ADDX (RB)+, R2 + RO = N*LN2~ -t0 + LN_FHI_LO 
50 «652.~=«C«60 0106 393 ADDX  R2, RO : RO=8 
B199 44) ; Compute A = N*LN2_HI + LN_FHI_HI and ALOG(x) | 
52 FEF2 CF 8E 65 0109 297 MULX3 (SP)+, LN2_HI, R2 : R2 = N®LN2_HI 
52 68 60 O10F 298 ADDX (RB), R2 ; R2 = A = N®LN2_HI + LN_FHI_HI 
50 52 60 0112 $33 ADDX R2, RO ; RO = A + B = ALOG(X) | 
50 0080 8F A2 0115 300 SUBW2  #1aX_EXP, RO : Divide by 2 
04 O11A 301 RET 
011B 0@ 
0118 39 NEG_EXP: 
56 55 A2 011B 304 SUBW sé, -R6 : RO = F = 2f 
56 56 QA OIE 305 MOVZBL R6, R6 3; R6 = index into ALOG table 
00000000 0121 306 IF NE, F_EXP=X_EXP 
6154 07 DIVW2 #1a<F_EXP=X_EXP>, R5 ; Shift RS to scale Bia 
0121 308 ENDC 
7E 55 6D O12) 09 CVTWX RS, =(SP) : Push N onto the st 
58 6846 «4098 «(0124 = 310 CVTBL (RB&)CR6], RB : RB = offset into FMI “sable | 
O1¢8 11 LN_1_PLUS WwW: 
50 19 0128 af BCSs LN_1_PLUS ; Branch to handle F close to 1 
58 00000000'GF48 7€£ O1gA 1 MOVAX  G*AMTAS$AB_D_FHICR8], R8 ; RB = Address of FHI | 
56 68 7D 01 é 1 MOVX (RB), R6 : R6 = FHI 
0138 16 : Compute Z = (F = FHI)/(F + FHI) 
135 17: 2 = EC Tex) = ¥"(V=X)I/0014X) + ¥"(1-xK)] 
Oise 18 : Z=(C1+X=- ¥ + xevJ/[1 + X + ¥ = X*¥) 
4 : }? : where Y = FHI®*2**N, roughly equal to (1+#X)/(1-x) 
56 55 «OA 61 .. 2a ADDW = RS, -RO : R6 = FHI * aeen = SFHI 
52 08 56 63 0138 é SUBX3 R6, $*#1.0, R2 : R2 = 1 = SFH 
56 08 60 013¢ ADDX $481.0, R6 : R6 = 1 + SFH 
3 56 66) (013F 4 DIVX R6, R : R2 = (1-SFHI)/(1+SFHI) = D 
56 50 «61 «0142 5 ADDX3 RO. Ro, R6 >: R6=D +X 
5 50 64 (0146 6 MULX RO, R > R2=D* X 
52 08 60 0149 7 ADDX $481.0, R2 : R02 = 1 + Dex 
56 52066 C4 5 DIVX  R2, R6 > RE = (D+X)/(14D*X) = 2 
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14F 9; 
be 9 ; Compute Z**2, P(Z**2) and Z*P(7**2) 
50 56 56 14F ; MULX3 R6, R6, RO ; RO = 2**2 | 
FEFF CF 05 29 $2 1 3 g POLYX Re. #LOGLENZ, LOGTAB2 ; RO = P(2*#2) 
50 6 64 ! : : MULX R6, R 3 RO = Z*P(Z**2) 
15¢ $ ; Compute B = N*LN2_LO + LN_FHI_LO + Z*P(Z*2) 
52 FEA7 CF 6E€ 65 01 Be MULKX3 (SP), LN2_LO, R2 ; R2 = N*LN2_LO 
52 7 60 016 9 ADDX =(RBS NS : RS = N*LN2-LO + LN_FHI_LO 
50 «652.~—=CG g16 40 ADDX  R2, RO > RO = 8 
O18 4g ; Compute A = N*LN2_HI + LN_FHI_HI and ALOG(X) 
52 FE93 CF BE 65 bre 44 MULX3 (SP)+#, LN2_HI, R2 ; R2 = N*LN2_HI 
26 78 62 O16€ 45 SUBX -(RB), R2 ; R2 = A = N®LN2_HI + LN_FHI_HI 
0 $2 60 0171 46 ADDX R2, R ; RO = A + B = ALOG(X) 
50 0080 BF A2 0174 47 SUBW2 #1aX_EXP, RO 3 Divide by 2 
04 0179 48 RET 
O17A 49 
O17A 50 ; 
O17A 51 ; Special logic for F close to 1 
O17A 26 ; 
O17A 5 
017A 54 LN_1_PLUS: 
56 08 50 63 OI7A 55 SUBX3 RO, S*#1.0, R6 3; R6 = 1-x 
55 B65 OI7E 56 TSTW R5 ; Determine which way to calculate W | 
OF 13 0180 57 BEQL fas 
56 10 56 67 0182 58 DIVX3 R6 sae. 0. R6 3 R6 = 2/(1-Xx) | 
56 08 62 0186 59 SUBX  $*#1 30 ‘ t RO = (14K)/ (1x) | 
56 55) a2 sO0189 Ss 360 SUBW SséRRS, : Scale R6 
56 08 62 018C 61 SUBX S*#1. 7 R6 ; RO =W 
09 11 O18F 6g BRB 20$ 
56 50 56 67 0191 63 10$: DIVX3 R6 at _ s R6 = / (1=xX) 
56 0080 8F AO 0195 364 ADDW = «#18x R6 i Re = Us Sex C1=X) = (14XD/CEXD = 1 | 
FE70 CF 98 56 75 O19A 65 208: POLYX R6, BLOGLENT LOGTAB1 ; RO = aww 
0 56 64 BiA0 66 MUL X Ro, RO : Finish B conput ing Q(wW) 
50 56 64 OIA 67 MULX RO ; RO = | 
52. FESD CF $f 65 0106 368 MULX3 86) LN2_LO, R2 > R2 = NSLN N2_LO 
50 52 60 O1AC 369 ADDX  R2, RO > RO = N®*LNO7LO + WeQ(W) 
50 56 60 OJ1AF 70 ADDX R6, RO ; RO = N*LN ay + wWeQ(W) + W 
6£ EEGA CF 64 01B¢ 71 MULX LNO_HI, (SP) : (SP) = NLR 
5 ge 6 18 if ADDX SP)+, R ; RO = ALOG (x } 
50 0080 BF A2 O1BA 37 SUBW2 #lax_EXP, RO : Divide by 2 
04 O1BF 74 RET 
1¢0 75 ; - 
188 6 ; X <= 0.0, signal error 
7E OO'8F 9A 1¢8 78 ERROR: RovzeL te INVARGMAT, -(SP) ; condition value 
50 01 OF 79 OQiC4 79 ASHQ #15, #T, RO ; RO = result = reserved operand -0.0 
1¢ 0 : goes © signal mechanism vector 
1¢ 1 s Crest’ MCH_RO/R1) so error handler 
1C ¢ 3 con modify the result. 
00000000'GF 01 FB OIC CALLS #1, G*MTHSSSIGNAL 3 tignel error and use real user's PC 
1CF 4 : ependent of CALL vs J 
04 OCF 5 RET : brat a - RO restored from 
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MTHSDATANH 
Symbol table 6-SEP=1984 11:21: 


The working set Limit was 1050 popes é k 

6942 bytes (14 pages) of virtual memory were used to buffer the intermediate code. 

There were 10 pages of symbol table space allocated to hold 20 non-local and 2 local symbols. 
447 source Lines were read in Pass 1, producing 11 object records in Pass 

3 pages of virtual memory were used to define 2 macros. 


AX/VMS Macro Vv04-00 


, g 
; Double Precision Hyperbolic Arctangent 16-SEP-1984 1:13:06 alyes os Bs 9 OP 


RD NOWRT NOVEC BYTE 
RD NOWRT NOVEC LONG 


ACMASK = 000001F 
ERR BS 0088 R 4 
ERROR 49 2 1 
iJ = 00 83 
LN2_H1 4464 4 & 01 
LN2_LO 8) 00 : R 01 
LN_T_PLUS 444d AR 01 
LN 1 PLUS_W 00 a3 R 01 
OGLEN1 = 000000 
LOGLEN2 = 0000000 
TABI 00000010 R 01 
LOGTAB2 00000058 R 4 
MTHSSAB_ALOG V teenreere xX 0 
MTHSS FHT eeereree XX 00 
MTHSSSIGNAL rerreere x 00 
MTHSDATANH 00000088 RG 01 
MTHSK_INVARGMAT eeeeeeee x 00 
NEG_ERP 0000011B R 01 
= 00000004 
X_EXP = 00000007 
pore wmoe ester en eee se + 
: Psect synopsis ! 
PSECT name Allocation PSECT No. Attributes 
es, 00000000 <( 0.) 00 ¢ OO.) NOPIC USR CON ABS ~ LCL NOSHR NO 
_MTHSCODE 000001D0 (¢ 464.) O1¢ 1.) PIC USR CON REL LCL 
gene een eee ee eer two mere e 
! Performance indicators ! 
ch ela cicarvcchiint-Iah atnc ata | 
Phase Page faults CPU Time Elapsed Time 
Initialization 30 ‘00:00:00.09 00:00:00.73 
Command processing 157 BO: B92 00 £6 00:00:04.94 
ass 101 00:00:01.3 00:00:05.19 
Symbol table sort 0 00:00:00.01 Be 20208 -01 
Pass 2 8s 00:00:00.93 00:00:03.49 
Symbol table output 00:00:00.04 00:00:00.12 
Psect synopsis output 2 00:00:00.02 00:00:00.02 
Cross-reference output 9 88: 8888-98 88: 88:00 -9? 
Assembler run totals 37 0:00:03. 0:00:14.54 
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toermeee renee wen ren oe eae weed 


! Macro Library Statistics ! 


jeuseneneneseerem eeeenr eres 


Macro Library name Macros defined 
_$255$DUA28: (SYSLIBISTARLET.MLB; 2 ante 8-9 98 
0 GETS were required to define 0 macros. 

There were no errors, warnings or information messages. 


MACRO/ENABLE=SUPPRESSION/DI SABLE=(GLOBAL , TRACEBACK) /LIS=LIS$:MTHDATANH/OBJ=OBJ$:MTHDATANH MSRC$:MTHJACKET/UPDATE=(ENHS:MTHJACKET) #MS 
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